# -*- coding: utf-8 -*-  
'''
PCA降维，主成分分析
    注：仅仅是练手。并没有对大数据做batch处理
    对于小数据量可以做到秒出结果，但大数据量算不出来。。。

Created on 2021年9月16日
@author: luoyi
'''
import numpy as np


#    PCA降维相关
class PCA:
    '''PCA降维：
        目标：让降维后的数据具有以下性质：
            1 同维度下方差最大（尽量保持数据差异）
            2 不同维度下协方差最小（尽量减小维度间依赖）
            用协方差矩阵来描述，其最好是个对角阵
        推导过程：
            假设：X为待降维的数据集，M维。先对X做标准化，每个维度均值为0：X - μ（方便后面求协方差矩阵）
                 Y为降维后的数据，K维
                     Y由X通过线性变换得到：Y = X * W    
                     W为变换矩阵：[M, K]
            问题即演变为：求W，使得Y的协方差矩阵对角线元素最大，非对角线元素为0
            设：
                X的协方差矩阵为∑x：
                    ∑x[i,j] = ∑(X)(xi - μi)(xj - μj) / N-1
                            = ∑(X)xi * xj / N-1                        ------ X事先标准化，均值为0
                Y的均值向量μ = { ∑(X)yi / N }
                           = { ∑(X)xi*W[i] / N }                       ------ W[i]表示W的第i列
                           = { W[i] * ∑(X)xi / N }
                           = { W[i] * μi }                             ------ X事先标准化，均值为0
                           = { 0 }
                Y的协方差矩阵为∑y：
                    ∑y[i,j] = ∑(Y)(yi - μi)(yj - μj) / N-1
                            = ∑(Y)(yi * yj) / N-1
                            = ∑(X)(xi*W[i] * xj*W[j]) / N-1
                            = ∑(X)W[i] * (xi * xj) * W[j] / N-1
                            = W[i] * ∑(X)(xi * xj)/N-1 * W[j]
                            = W[i] * ∑x[i,j] * W[j]
                         ∑y = W.T * ∑x * W
            则：目标为找到矩阵W，让∑y变为对角化。条件是∑x是实对称阵。。。
            眼熟不，形似P.T * A * P = V
            其中：
                A为实对称阵
                P为A的特征向量组成的矩阵，且P的每行两两正交
                V为每个P对应的特征值
            那么W就取为∑x的特征值组成的矩阵[M,M]，且按特征值从大到小排序
            舍弃最小的M-K个特征值得到的矩阵W[M,K]，即为所求变换矩阵。其满足性质：
                1 ∑y是对角阵
                2 ∑y的主对角线元素保留下来的都是最大的前K个特征值（在保证条件1的基础上尽量满足条件2）
    '''
    def __init__(self, K):
        '''
            @param K: 将第到K维
        '''
        self._K = K
        pass
    
    #    降维
    def reduce(self, X):
        '''
            @param X: Ndarray(N, M)    需要降维的数据。N：数据长度，M：数据维度
            @return: Y Ndarray(N, K)    降维后的数据
        '''
        #    均值变为0
        X = self.mean_zero(X)
        
        #    X的协方差矩阵
        cov_x = np.cov(X.T)
        
        #    协方差矩阵求特征值和对应的特征向量
        _, p = np.linalg.eig(cov_x)
        
#            取前K个特征向量[M, K]
        p = p[:, :self._K]
        
        #    计算降维后的Y
        Y = np.dot(X, p)
        return Y
    #    标准化
    def mean_zero(self, X):
        mean = np.mean(X, axis=0)
        return X - mean
    pass



